#include "analyzer_common.h"

void exclusive_scan_sparse_pointer(sparse_pointer_t *input, int length)
{
    if (length == 0 || length == 1)
        return;

    sparse_pointer_t old_val, new_val;

    old_val = input[0];
    input[0] = 0;
    for (int i = 1; i < length; i++)
    {
        new_val = input[i];
        input[i] = old_val + input[i - 1];
        old_val = new_val;
    }
}

int findlevel(
    const sparse_pointer_t *cscColPtr,
    const sparse_index_t *cscRowIdx,
    const sparse_pointer_t *csrRowPtr,
    const sparse_index_t m,
    int *nlevel,
    int *levelPtr,
    int *levelItem)
{
    // prepare arrays for level-sets of size maximum m or n
    int *indegree = (int *)malloc(m * sizeof(int));

    // prepare in-degree
    for (int i = 0; i < m; i++)
    {
        indegree[i] = csrRowPtr[i + 1] - csrRowPtr[i];
    }
    // for (int i = 0; i < m; i++)
    //     printf("%d ", indegree[i]);
    // for (int i = 0; i < m; i++)
    // printf("%d ", indegree[i]);
    // printf("\n");


    // find root items
    int lv = 0;
    int ptr = 0;

    levelPtr[0] = 0;

    for (int i = 0; i < m; i++)
    {
        if (indegree[i] == 1)
        {
            levelItem[ptr] = i;
            ptr++;
        }
    }

    // for (int i = 0; i < ptr; i++)
    // printf("%d ", levelItem[i]);
    // printf("\n");
    
    // #items in the 1st level
    levelPtr[1] = ptr;
    // printf("shoule lvptr = %d\n\n\n", levelPtr[1]);

    int lvi = 1;
    while (levelPtr[lvi] != m)
    {
        for (int i = levelPtr[lvi - 1]; i < levelPtr[lvi]; i++)
        {
            int node = levelItem[i];
            for (int j = cscColPtr[node]; j < cscColPtr[node + 1]; j++)
            {
                int visit_node = cscRowIdx[j];
                indegree[visit_node]--;
                if (indegree[visit_node] == 1)
                {
                    levelItem[ptr] = visit_node;
                    ptr++;
                }
            }
        }
        lvi++;
        levelPtr[lvi] = ptr;
    }

    *nlevel = lvi;
        // printf("lvi = %d\n", lvi);
    // printf("levelPtr:\n");
    // for (int i = 0; i < m; i++)
    //     printf("%d ", levelPtr[i]);

    free(indegree);

    return 0;
}

void matrix_transposition(
    const sparse_index_t m,
    const sparse_index_t n,
    const sparse_pointer_t nnz,
    const sparse_pointer_t *csrRowPtr,
    const sparse_index_t *csrColIdx,
    const sparse_value_t *csrVal,
    sparse_index_t *cscRowIdx,
    sparse_pointer_t *cscColPtr,
    sparse_value_t *cscVal)
{
    // histogram in column pointer
    memset(cscColPtr, 0, sizeof(sparse_pointer_t) * (n + 1));
    for (int i = 0; i < nnz; i++)
    {
        cscColPtr[csrColIdx[i]]++;
    }

    // prefix-sum scan to get the column pointer
    exclusive_scan_sparse_pointer(cscColPtr, n + 1);

    // for (int i = 0; i < m + 1; i++)
    //     printf("%d ", cscColPtr[i]);
    // printf("\n\n");

    sparse_pointer_t *cscColIncr = (sparse_pointer_t *)malloc(sizeof(sparse_pointer_t) * (n + 1));
    memcpy(cscColIncr, cscColPtr, sizeof(sparse_pointer_t) * (n + 1));

    // for (int i = 0; i < n + 1; i++)
    //     printf("%d ", cscColIncr[i]);
    // printf("\n\n");

    // insert nnz to csc
    for (int row = 0; row < m; row++)
    {
        for (int j = csrRowPtr[row]; j < csrRowPtr[row + 1]; j++)
        {
            int col = csrColIdx[j];

            cscRowIdx[cscColIncr[col]] = row;
            // printf("%d ", cscRowIdx[cscColIncr[col]]);
            cscVal[cscColIncr[col]] = csrVal[j];
            cscColIncr[col]++;
        }
        // printf("\n\n");
    }

    free(cscColIncr);
}

int findlevel_csr(const sparse_pointer_t *csrRowPtr,
                  const sparse_index_t *csrColIdx,
                  const sparse_value_t *csrVal,
                  const sparse_index_t m,
                  const sparse_index_t n,
                  const sparse_pointer_t nnz,
                  int *nlevel,
                  int *parallelism_min,
                  int *parallelism_avg,
                  int *parallelism_max)
{
    if (m != n)
    {
        printf("This is not a square matrix, return.\n");
        return -1;
    }
    struct timeval start_time;

    // transpose to have csr data
    sparse_pointer_t *cscColPtr = (sparse_pointer_t *)malloc((n + 1) * sizeof(sparse_pointer_t));
    sparse_index_t *cscRowIdx = (sparse_index_t *)malloc(nnz * sizeof(sparse_index_t));
    sparse_value_t *cscVal = (sparse_value_t *)malloc(nnz * sizeof(sparse_value_t));

    analyzer_timer_start(&start_time);
    printf("[findlevel] start transposition\n");
    // transpose from csc to csr
    matrix_transposition(m, n, nnz,
                         csrRowPtr, csrColIdx, csrVal,
                         cscRowIdx, cscColPtr, cscVal);
    printf("[findlevel] end transposition %.2lf s\n", analyzer_timer_stop(&start_time));

    int *levelPtr = (int *)malloc((m + 1) * sizeof(int));
    int *levelItem = (int *)malloc(m * sizeof(int));

    analyzer_timer_start(&start_time);
    printf("[findlevel] start findlevel\n");
    int nlv = 0;
    findlevel(cscColPtr, cscRowIdx, csrRowPtr, m, &nlv, levelPtr, levelItem);
    printf("[findlevel] end findlevel %.2lf s\n", analyzer_timer_stop(&start_time));

    // calculate min, avg, and max #item in levels
    int min = levelPtr[1] - levelPtr[0];
    int max = levelPtr[1] - levelPtr[0];
    for (int i = 1; i < nlv; i++)
    {
        int nitem = levelPtr[i + 1] - levelPtr[i];
        if (nitem > max)
            max = nitem;
        else if (nitem < min)
            min = nitem;
    }
    *parallelism_min = min;
    *parallelism_max = max;
    *parallelism_avg = m / nlv;

    *nlevel = nlv;

    free(levelPtr);
    free(levelItem);
    free(cscColPtr);
    free(cscRowIdx);
    free(cscVal);

    return 0;
}

int level_set_detect(sparse_matrix_t* s){

    int nlevel, parallelism_min, parallelism_max, parallelism_avg;

    struct timeval start_time;
    analyzer_timer_start(&start_time);
    findlevel_csr(
        s->rowpointer, 
        s->columnindex,
        s->value,
        s->m,
        s->n,
        s->nnz,
        &nlevel,
        &parallelism_min,
        &parallelism_avg,
        &parallelism_max
    );
    printf("level max = %d, min = %d, avg = %d nlevel = %d %.2lf s\n", parallelism_max, parallelism_min, parallelism_avg, nlevel, analyzer_timer_stop(&start_time));
}